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Abstract 

We consider a clustered network with small-world sub-networks of inhibitory fast spiking in¬ 
terneurons, and investigate the effect of inter-modular connection on emergence of fast sparsely 
synchronized rhythms by varying both the inter-modular coupling strength Jinter and the average 
number of inter-modular links per interneuron ■ In contrast to the case of non-clustered 

networks, two kinds of sparsely synchronized states such as modular and global synchronization 
are found. For the case of modular sparse synchronization, the population behavior reveals the 
modular structure, because the intra-modular dynamics of sub-networks make some mismatch¬ 
ing. On the other hand, in the case of global sparse synchronization, the population behavior is 
globally identical, independently of the cluster structure, because the intra-modular dynamics of 
sub-networks make perfect matching. We introduce a realistic cross-correlation modularity mea¬ 
sure, representing the matching-degree between the instantaneous sub-population spike rates of the 
sub-networks, and examine whether the sparse synchronization is global or modular. Depending on 
its magnitude, the inter-modular coupling strength Jinter seems to play “dual” roles for the pacing 
between spikes in each sub-network. For large Jinter-, due to strong inhibition it plays a destructive 
role to “spoil” the pacing between spikes, while for small Jinter it plays a constructive role to “fa¬ 
vor” the pacing between spikes. Through competition between the constructive and the destructive 
roles of Jinter, there exists an intermediate optimal Jinter at which the pacing degree between spikes 
becomes maximal. In contrast, the average number of inter-modular links per interneuron 
seems to play a role just to favor the pacing between spikes. With increasing the pacing 

degree between spikes increases monotonically thanks to the increase in the degree of effectiveness 
of global communication between spikes. Furthermore, we employ the realistic sub- and whole- 
population order parameters, based on the instantaneous sub- and whole-population spike rates, 
to determine the threshold values for the synchronization-unsynchronization transition in the sub- 
and whole-populations, and the degrees of global and modular sparse synchronization are also mea¬ 
sured in terms of the realistic sub- and whole-population statistical-mechanical spiking measures 
defined by considering both the occupation and the pacing degrees of spikes. It is expected that 
our results could have implications for the role of the brain plasticity in some functional behaviors 
associated with population synchronization. 

PACS numbers: 87.19.1m, 87.19.1c 
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I. INTRODUCTION 


Recently, much attention has been paid to brain rhythms in health and disease mi. 
Particularly, we are interested in fast sparsely synchronized cortical rhythms which are as¬ 
sociated with diverse cognitive functions such as sensory perception, feature integration, 
selective attention, and memory formation [3]. At the population level, local held potential 
recordings have been observed to show synchronous fast oscillations [e.g., gamma rhythm 
(30-100 Hz) and ultrafast sharp-wave ripple (100-200 Hz)], while individual neuronal record¬ 
ings have been found to exhibit stochastic and intermittent spike discharges [iflTU]. Thus, 
single-cell bring activity differs markedly from the population oscillatory behavior. These 
sparsely synchronized rhythms are in contrast to fully synchronized rhythms. For the case 
of full synchronization, individual neurons hre regularly at the population frequency like the 
clock oscillators m- Hence, the fully synchronized oscillations may be well described by 
using the conventional coupled-oscillator model composed of suprathreshold spiking neurons 
above a threshold in the absence of noise or for weak noise [12]. However, such coupled- 
oscillator models are not adequate for describing sparse synchronization because individual 
neurons hre stochastically at low rates like the Geiger counters. Brunei et ah in [T31 - 1T8] 
developed a framework appropriate for description of fast sparse synchronization by taking 
an opposite view from that of coupled oscillators. Under the condition of strong external 
noise, suprathreshold spiking neurons discharge irregular brings as Geiger counters, and then 
the population state becomes unsynchronized. However, when inhibitory recurrent feedback 
becomes sufficiently strong, this asynchronous state may be destabilized, and then a syn¬ 
chronous population state with irregular and intermittent individual discharges emerges. 
For this case, average total (external excitatory plus recurrent inhibitory) input current into 
individual neurons is subthreshold, but stochastic and intermittent brings are triggered when 
buctuations (due to noise in external and recurrent inputs) cross a threshold. In this way, 
under the balance between strong external noise and strong recurrent inhibition, fast sparse 
synchronization was found to occur in networks of suprathreshold neurons HMEj. Similar 
sparsely synchronized rhythms were also found to appear via cooperation of noise-induced 
spikings of subthreshold neurons (which can not bre spontaneously without noise) [T91 - I2T] . 
However, in contrast to the above works on suprathreshold neurons, sparse synchronization 
for the case of subthreshold neurons has been found to appear under relatively weak external 
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noise and recurrent inhibition, and the sparsely synchronized rhythms were also found to be 
slow when compared with the suprathreshold case. 

In this paper, we are concerned about emergence of fast sparsely synchronized rhythms 
in an ensemble of suprathreshold neurons, as in the previous works of Brunei et ah [IMH] 
where both random and global synaptic couplings were considered. However, connection 
architecture of the real brain has been found to have complex topology which is neither reg¬ 
ular nor random [221 - I30] . Particularly, mammalian (e.g., cat and macaque) brain anatomical 
networks and human brain functional (fMRI) networks have been revealed to have a mod¬ 
ular structure composed of relatively sparsely linked clusters with spatial localization, as in 
social, technological, and biological complex networks [MI - 138] . Within each cluster, nodes 
are highly inter-connected and exhibit similar connectional and functional features. This 
clustered organization of the brain network reveals the anatomical substrate for segregation 
which refers to the subdivision of the brain into regions specialized in particular functional 
tasks [MHH]- This segregation allows the brain to process information in parallel, simulta¬ 
neously by distinct populations of neurons. However, for emergence of a coherent perception 
and comprehensive understanding of the environment as a whole, specialized information 
of different modalities and features should be integrated. This integration refers to the ca¬ 
pacity of a system to collect information of different nature and combine it to produce new 
useful information. For example, sensory perception requires the binding of the features 
of a receptive field (e.g., color, orientation, and position of a visual object). In this way, 
brain connectivity should be organized into a balance between segregation (specialization) 
and integration (binding) jlSHlG]. Here, in our brain network we take into consideration 
the modular structure of the real brain. For modeling the modular structure of real brain, 
we consider a clustered network composed of sparsely connected sub-networks. The sub¬ 
networks stand for the modules (clusters) of the brain. As is also known, the connection 
structure in each module of the real brain reveals complex topology such as small-worldness 
and scale-freeness [221130] . Here, each sub-network (representing a cluster) is modeled as 
the Watts-Strogatz small-world network which interpolates between the regular lattice with 
high clustering (i.e., high cliquishness of a typical neighborhood) and the random graph 
with short path length (i.e., average short separation between two neurons represented by 
average number of synapses between two neurons along the minimal path) by varying the 
rewiring probability Prewiring from local to long-range connections; Prewiring = 0 and 1 corre- 
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spond to the regular lattice and the random graph, respectively |13-ll9] . The Watts-Strogatz 
model for the small-world sub-network may be regarded as a cluster-friendly extension of 
the random network by reconciling the six degrees of separation (small-worldness) [SDl El] 
with the circle of friends (clustering). Many recent works on various subjects of neurody¬ 
namics have been done in small-world networks with predominantly local connections and 
rare long-distance connections [52fl63] . Effect of this small-world connectivity on fast sparse 
synchronization has also been studied in our recent work [M] . 

We note that real brain networks, consisting of sparsely inter-connected modules, are 
far more complex than minimal non-modular models such as small-world and scale-free 
networks. The main purpose of our study is to investigate emergence of sparsely synchro¬ 
nized rhythms in more realistic modular networks. Clustered neural networks composed of 
regular, small-world, and scale-free sub-networks have been employed for study on several 
subjects of neurodynamics [391 - 1111 I65l - t69] . For our aim, we consider a clustered network 
with small-world sub-networks of inhibitory spiking neurons, and investigate the effect of 
inter-modular connection on emergence of fast sparsely synchronized rhythms. In the ab¬ 
sence of inter-modular coupling, we consider three cases for the intra-modular dynamics 
in sub-networks: (1) synchronized in all identical sub-networks, (2) unsynchronized in all 
identical sub-networks, and (3) synchronized/unsynchronized in non-identical sub-networks. 
For each case, we study the population states by changing both the inter-modular coupling 
strength Jinter and the average number of inter-modular links per interneuron Con¬ 

sequently, two kinds of sparse synchronization such as modular and global synchronization 
are found, in contrast to the case of non-modular networks HMHl El EO]. For the case of 
modular sparse synchronization, the population behavior reveals the modular structure, and 
hence the degree of sparse synchronization in the whole population becomes less than that 
in the sub-networks, because the intra-modular dynamics of sub-networks make some mis¬ 
matching. In contrast, for the case of global sparse synchronization, the population behavior 
is globally identical, independently of the cluster structure, because the intra-modular dy¬ 
namics of sub-networks make perfect matching. These modular and global synchronization 
may be well visualized in the raster plots of spikes. For the case of synchronization, syn¬ 
chronous “stripes” (composed of spikes and indicating population synchronization) appear 
successively in the raster plots, while spikes are completely scattered (without forming any 
stripes) in the case of unsynchronization. Synchronization pacing (representing the smear- 
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ing of spiking stripes: less smearing, better pacing) varies depending on Jinter- For large 
■Jinter F pl^YS a destructive role to spoil the pacing between sparse spikes, because of strong 
inhibition. Hence, when passing a large threshold a transition to unsynchronization occurs. 
However, for small Jinter it plays a constructive role to favor the pacing between spikes in 
each sub-network. Hence, via competition between the constructive and the destructive roles 
of Jinter) there appears an intermediate optimal Jinter at which the pacing degree between 
spikes becomes maximal. In this way, Jinter plays dual roles for the pacing between spikes 
in each sub-network, depending on its magnitude. On the other hand, the average number 
of inter-modular links per interneuron plays a role just to favor the pacing between 

spikes. As is increased, the pacing degree between spikes increases monotonically 

due to the increase in the degree of effectiveness of global communication between spikes. 
To make characterization of the synchronization-unsynchronization transitions in the sub- 
and whole-populations, we employ the realistic sub- and whole-population order param¬ 
eters, based on the instantaneous sub- and whole-population spike rates HU. Moreover, 
we introduce a realistic cross-correlation modularity measure, representing the matching- 
degree between the instantaneous sub-population spike rates of sub-networks, and examine 
whether the sparse synchronization is global or modular. The degrees of modular and global 
sparse synchronization are also measured in terms of the realistic sub- and whole-population 
statistical-mechanical spiking measures dehned by considering both the occupation and the 
pacing degrees of the spikes mi. 

This paper is organized as follows. In Sec. we describe a clustered network with small- 
world sub-networks of inhibitory fast spiking (FS) interneurons, and then the governing 
equations for the population dynamics are given. Detailed explanations on methods for 
characterization of individual and population states in clustered networks are also given in 
Sec. |III[ Then, in Sec. W we investigate the effect of inter-modular connection on emergence 
of fast sparsely synchronized rhythms by varying both Jinter and Finally, a summary 

is given in Section |V| 
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II. CLUSTERED SMALL-WORLD NETWORK OF INHIBITORY FS IZHIKE- 


VICH INTERNEURONS 


In this section, we first describe our cluster network composed of M(= 3) small-world 


sub-networks, each of which contains L FS Izhikevich interneurons in the subsection IIA 


Then, the governing equations for the population dynamics in the clustered small-world 


network are given in the subsection IIB 


A. Clustered Small-World Networks 

We consider a clustered network with M (= 3) small-world sub-networks. Each small- 
world sub-network consists of L inhibitory interneurons equidistantly placed on a one¬ 
dimensional ring of radius L{27r. For illustrative purpose, an example of the clustered 
network topology is shown in Fig. [Tj Each of the three sub-networks, consisting of L = 20 
interneurons, is modeled as the Watts-Strogatz small-world network which interpolates be¬ 
tween the regular lattice and the random graph by varying the rewiring probability Prewiring 
from local to long-range connections We start from the case of Prewiring = 0, cor¬ 

responding to a directed regular ring lattice where each interneuron is coupled to its hrst 
(= 4) neighbors /2 on either side) via outward synapses. Then, we rewire 

each outward connection at random with probability Prewiring such that self-connections and 
duplicate connections are excluded, and the value of Prewiring is 0.25 for the case of Fig. 
Within each small-world sub-network, the average number of intra-modular synaptic inputs 
per interneuron is while there exist 8 sparse random inter-modular links between 

small-world sub-networks. 


B. Governing Equations for The Population Dynamics 

As an element in our clustered small-world network, we choose the FS Izhikevich interneu¬ 
ron model which is not only biologically plausible, but also computationally efficient HMg. 
We consider the clustered network composed of M(= 3) small-world subnetworks, each of 
which consists of L FS interneurons; L = 10^, except for the case of order parameters and 
spatial cross-correlation functions. The following equations Q-Q govern the population 




dynamics in the clustered small-world network: 
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Here, vf\t) and u'l\t) are the state variables of the ith interneuron in the Jth sub-network 
at a time t which represent the membrane potential and the recovery current, respectively. 
These membrane potential and the recovery variable, vl^\t) and u\^\t), are reset according 
to Eq. (3) when vl^\t) reaches its cutoff value Vp. C, Vr, and Vt in Eq. (j^ are the membrane 
capacitance, the resting membrane potential, and the instantaneous threshold potential, 
respectively. The parameter values used in our computations are listed in Table [T| More 
details on the FS Izhikevich interneuron model, the external stimulus to the FS interneuron, 
the intra-modular and the inter-modular synaptic currents, and numerical integration of the 
governing equations are given in the following subsubsections. 


1. FS Izhikevich Interneuron Model 

The Izhikevich model matches neuronal dynamics by tuning the parameters {k, a, 6, c, d) 
instead of matching neuronal electrophysiology, in contrast to Hodgkin-Huxley-type 
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conductance-based models. The parameters k and b are related to the neuron’s rheobase 
and input resistance, and a, c, and d are the recovery time constant, the after-spike reset 
value of V, and the total amount of outward minus inward currents during the spike and 
affecting the after-spike behavior (i.e., after-spike jump value of u), respectively. Depend¬ 
ing on the values of these parameters, the Izhikevich neuron model may exhibit 20 of the 
most prominent neuro-computational features of cortical neurons [72HZ5]. Here, we use the 
parameter values for the FS interneurons in the layer 5 rat visual cortex, which are listed in 
the 1st item of Table U 

2. External Stimulus to The FS Izhikevich Interneuron 

Each Izhikevich interneuron is stimulated by both a common DC current I^c and an 
independent Gaussian white noise as shown in the 3rd and the 4th terms in Eq. (j^. 
The Gaussian white noise satishes = 0 and = 6ij Sij — where 

(■ ■ ■) denotes the ensemble average. Here the Gaussian noise is a parametric one that 
randomly perturbs the strength of the applied current Idc, and its intensity is controlled by 
the parameter D. For D = 0, the Izhikevich interneuron exhibits a jump from a resting state 
to a spiking state via subcritical Hopf bifurcation for a higher threshold Idc, hi— 73.7) by 
absorbing an unstable limit cycle born via a fold limit cycle bifurcation at a lower threshold 
Idc,i{— 72.8). Therefore, the Izhikevich interneuron shows type-H excitability since it begins 
to hre with a non-zero frequency [76l[77]. With increasing Idc from lDC,h, the mean firing 
rate / increases monotonically. The values of Idc and D used in this paper are given in the 
2nd item of Table HI 


3. Intra-modular and Inter-modular Synaptic Currents 


The last two terms in Eq. Q represent the intra- and the inter-modular synaptic cou¬ 
plings of inhibitory FS interneurons, jdnter,syn) represent 

the intra- and the inter-modular synaptic currents injected into the Hh neuron in the Jth sub¬ 
network, respectively. The synaptic connectivity is given by the connection weight matrix 
W i={w^j’'^'^}) where = 1 if the neuron j in the Jth sub-network is presynaptic to the 

neuron i in the Jth sub-network; otherwise, = 0. Here, the intra-modular synaptic con- 
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nection is modeled in terms of the Watts-Strogatz small-world network. Then, the in-degree 
of the ith neuron in the Jth sub-network for the intra-modular synaptic connection, 

(i.e., the number of intra-modular synaptic inputs to the neuron i in the Jth sub-network) is 
given by For this intra-modular case, the average number of intra- 

modular synaptic inputs per neuron is 'I2i=i In contrast to the 

intra-modular connection, the inter-modular synaptic connection is given randomly. Then, 
the in-degree of the zth neuron in the Jth sub-network for the inter-modular synaptic con¬ 
nection, dyy*®'’ (i.e., the number of inter-modular synaptic inputs to the neuron i in the Jth 
sub-network) is given by In the inter-modular case, the aver¬ 
age number of inter-modular synaptic inputs per neuron is Yld=\ Yl!i=i 

and these inter-modular links are randomly connected with the inter-modular connection 

j^(inter) 

probability Pinter = Compared to the intra-modular connections, the inter-modular 

connections are sparse (i.e., The fraction of open synaptic ion chan¬ 

nels at time t is denoted by s{t). The time course of sf\t) of the jth neuron in the Jth 


sub-network is given by a sum of delayed double-exponential functions E(t — t 

Eq. (7)], where r/ is the synaptic delay, and and qj-q the /th spiking time and the 
total number of spikes of the jth neuron in the Jth sub-network at time t, respectively. Here, 
E{t) [which corresponds to contribution of a presynaptic spike occurring at time 0 to s{t) 
in the absence of synaptic delay] is controlled by the two synaptic time constants: synaptic 
rise time and decay time and 0(t) is the Heaviside step function: 0(f) = 1 for t > 0 
and 0 for f < 0. The intra- and the inter-modular synaptic coupling strengths are controlled 
by the parameters Jintra and Jinteri respectively and Vgyn is the synaptic reversal potential. 
For the inhibitory GABAergic synapse (involving the GABAa receptors), the values of r/, 
Tr, Td, and Vsyn are listed in the 3rd item of Table 


see 


4- Numerical Method 

Numerical integration of stochastic differential Eqs. 0-0 is done by employing the 
Heun method [TH] with the time step At = 0.01 ms. For each realization of the stochastic 
process, we choose a random initial point [^^^(0), ^^^(0)] for the ith (f = 1,..., 77) neuron 
in the Jth sub-network with uniform probability in the range of ^^^^^(0) G (—50, —45) and 
wf^o) e (10,15). 
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III. METHODS FOR CHARACTERIZATION OF INDIVIDUAL AND POPULA¬ 
TION STATES IN CLUSTERED NETWORKS 

In the following subsections, we explain methods used to characterize individual and pop¬ 
ulation states in clustered networks. Particularly, emergence of population synchronization 
and its degree are characterized by employing realistic measures, based on instantaneous 
sub- and whole-population spike rates HD. Furthermore, we introduce a realistic cross¬ 
correlation modularity measure, denoting the matching-degree between the instantaneous 
sub-population spike rates to examine whether the population synchronization is global or 
modular. 

A. Characterization of Individual Firing Behaviors 

Firing behaviors of individual interneurons are characterized in terms of the inter-spike 
interval (ISI) histogram and the mean bring rate (MFR) distribution. The ISI histogram 
is composed of 5 x 10^ ISIs (obtained from all the interneurons), and the bin size for the 
histogram is 0.5 ms. The MFR for each interneuron is calculated by following the membrane 
potential during the averaging time of 2 x 10"^ ms after discarding the transient time of 10^ 
ms, and the bin size for the histogram is 2 Hz. 

B. Sub- and Whole-population Variables 

In computational neuroscience, an ensemble-averaged sub-population potential 
for the Jth sub-network (/ = 1,2,3), containing L FS Izhikevich interneurons, 

F'’W = lET’W. (8) 

i=l 

and an ensemble-averaged whole-population potential for the whole network with M (=3) 
sub-networks, 

1 ^ 

^ I=l 

are often used for describing emergence of population neural synchronization in the sub- 
and the whole-populations, respectively (e.g., sparse synchronization in a population of 
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subthreshold neurons was described in terms of an ensemble-averaged global potential [19]- 
[21] )• However, to directly obtain Vs^^\t) and 14;(t) in real experiments is very difficult. To 
overcome this difficulty, instead of Vy\t) and 14 (4) experimentally-obtainable 

instantaneous sub- and whole-population spike rates which are often used as collective quan¬ 
tities showing sub- and whole-population behaviors [allinHIHlElllTollTI]. The instantaneous 
sub-population spike rate (ISPSR) Ri^\t) is obtained from the raster plot of neural spikes 
which is a collection of spike trains of individual neurons in the Jth sub-population. Such 
raster plots of spikes, where sub-population spike synchronization may be well visualized, 
are fundamental data in experimental neuroscience. For the synchronous case, “stripes” 
(composed of spikes and indicating sub-population synchronization) are found to be formed 
in the raster plot. Hence, for a synchronous case, an oscillating ISPSR appears, while for 
an unsynchronized case the ISPSR is nearly stationary. To obtain a smooth ISPSR, we 
employ the kernel density estimation (kernel smoother) [79|. Each spike in the raster plot is 
convoluted (or blurred) with a kernel function Kh{t) to obtain a smooth estimate of ISPSR 
for the Jth sub-network, Ri^\t): 

L 4 '^ 

= ( 10 ) 

i=l 

where is the sth spiking time of the Rh neuron in the Jth sub-network, nf'* is the total 
number of spikes for the Rh neuron in the Jth sub-network, and we use a Gaussian kernel 
function of band width h\ 


Khit) 



— oo < t < oo. 


( 11 ) 


Throughout the paper, the band width of the Gaussian kernel estimate is h = 1 ms. 
Then, the instantaneous whole-population spike rate (IWPSR) kernel estimate Rw{t) for 
the whole population is given by an average of the ISPSR kernel estimates of the M(= 3) 
sub-populations: 

1 ^ 

= ( 12 ) 

^ I=l 

Moreover, for the synchronous case, the sub-population frequency fp^'^ of the regularly- 
oscillating ISPSR Ri^\t) may be obtained from the one-sided power spectrum of ARi^\t) 
[= Ri^\t) — i?i'^^(f)] with the mean-squared amplitude normalization. The number of data 
for the power spectrum is 2^^, and the overline denotes the time average. 


13 




C. Sub- and Whole-population Order Parameters 


As is well known, a conventional order parameter, based on the ensemble-averaged global 
potential, is often used for describing transition from synchronization to unsynchronization 
in computational neuroscience [T9ll^ IMfl82] . Recently, instead of the global potential, we 
used an experimentally-obtainable instantaneous population spike rate kernel estimate, and 
developed a realistic order parameter for the case of the non-modular networks, which may 
be applicable in both the computational and the experimental neuroscience m- For the 
case of modular networks, the mean square deviation of the ISPSR kernel estimate E!'s\t) 
for the Jth sub-network (/=!, 2, 3), 

= (13) 

and the mean square deviation of the IWPSR kernel estimate Rw{t) for the whole network, 

( 14 ) 

play the role of realistic sub- and whole-population order parameters and to deter¬ 
mine a threshold for the synchronization-unsynchronization transition, where the overbar 
represents the time average. Here, each order parameter is obtained through average over 20 
realizations, and the averaging time for the calculation of the order parameter in each realiza¬ 
tion is 4 X 10^ ms. Then, the order parameters and representing the time-averaged 
fluctuations of Ri^\t) and Rw{t), approach non-zero (zero) limit values for the synchronized 
(unsynchronized) state in the thermodynamic limit of L —)■ oo. These order parameters 
may be regarded as thermodynamic measures because they concern just the macroscopic 
ISPSR and IWPSR kernel estimates Ri^\t) and Rw{t) without any consideration between 
the macroscopic ISPSR and IWPSR kernel estimates and microscopic individual spikes. 

D. Spatial Cross-correlation Functions 

To further understand the synchronization-unsynchronization transition, we consider the 
“microscopic” dynamical cross-correlations between neuronal pairs. For obtaining dynam¬ 
ical pair cross-correlations, each spike train of the ith neuron in the Jth sub-network is 
convoluted with a Gaussian kernel function Kh{t) of band width h to get a smooth estimate 
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of instantaneous individual spike rate (IISR) 

= (15) 

s=l 

where is the sth spiking time of the zth neuron in the Jth sub-network, is the total 
number of spikes for the ith neuron, and Kh{t) is given in Eq. 0 - Then, the normalized 
temporal cross-correlation function between the IISR kernel estimates rf\t) and 

r^p{t) of the (f, j) neuronal pair in the Jth sub-network is given by: 


ArP{t + T)ArP{t) 

A'‘Aw\/ 


( 16 ) 


where ArP{t) = rP(t) —rP{t) and the overline denotes the time average. Here, the number 
of data used for the calculation of each temporal cross-correlation function cP{t) is 4 x 10^. 
Similar to the case of non-modular small-world network [6l], we introduce the spatial cross¬ 
correlation function CP (/ = l,...,L/2) between neuronal pairs separated by a spatial 
distance I in the Jth sub-network through average of all the temporal cross-correlations 
between rP{t) and {i = 1, ...,L) at the zero-time lag: 

d'* = T E dk(0) for i = 1, ■ ■ ■ , r/2. (17) 

i=l 

Here, if i -|- / > L in Eq. 0, then i + l — L is considered instead of i -|- / because neurons lie 
on the ring. If the spatial cross-correlation function CP {I = 1, ...,L/2) is non-zero in the 
whole range of I, then the spatial correlation length rji becomes L/2 (note that the maximal 
distance between neurons is L/2 because of the ring architecture on which neurons exist) 
covering the whole sub-network. For this case, synchronization appears in the sub-network; 
otherwise, unsynchronization occurs. 


E. Cross-correlation Modularity Measure 

To determine the type of synchronization (modular or global) in modular networks, we 
measure the matching degree between the intra-modular dynamics of sub-networks in terms 
of the cross-correlation modularity measure Cm, based on ISPSRs. The normalized temporal 
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cross-correlation function Ci^j{t) between the ISPSR kernel estimates and of 

the Jth and the Jth sub-networks is given by: 


ARP{t + T)ARi^\t) 


(18) 


where ARi^\t) = Ri^\t) — Ri^\t) and the overline denotes the time average. Figures ^cl)- 
[^c5) show the normalized temporal cross-correlation functions C/j(r) for Rnter =10, 30, 70, 
400, and 1200, respectively. Then, the cross-correlation modularity measure Cm is obtained 
through average of the temporal cross-correlations between all the sub-population pairs at 
the zero-time lag: 


M-l M 


Cm = 


M(M- 1) 


rTi5Z C/.Xo). 


(19) 


7=1 J=7+l 


Here, the cross-correlation modularity measure {CM)r is obtained through average over 
20 realizations, and the number of data used for the calculation of each temporal cross¬ 
correlation function in each realization is 4 x 10^. 


F. State Diagram 

Population states vary depending on the inter-modular connection parameters Rnter and 
Ms^^\ which may be well shown in the state diagram in the Rnter — plane. To ob¬ 
tain the state diagram, we hrst divide the Rnter — plane into the 20 x 10 grids. Then, 

at each grid point, we calculate the sub-population order parameters (/ = 1, 2, 3) for 
L = 10^ and 10^ to determine whether the population state at the grid point is synchronized 
or unsynchronized. If for L = 10^ is smaller than / ■ for L = 10^ (/ is some appro¬ 
priate factor less than unity; for convenience we set / = 0.3), is expected to decrease 
with increasing L. For the case of decrease in with increasing L, unsynchronization 
occurs at the grid point; otherwise, synchronization appears. Next, at the grid points where 
synchronization occurs, we calculate the cross-correlation modularity measure Cm to deter¬ 
mine whether the population synchronization is modular or global. If Cm is larger than a 
threshold (in our computation, we set Cj^} = 0.995), global synchronization is expected 
to appear at the grid point; otherwise, modular synchronization occurs. After determin¬ 
ing the population states (modular or global synchronization and unsynchronization) at all 
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grid points, we try to obtain the synchronization-unsynchronization and the modular-global 
synchronization transition curves accurately. To this end, we calculate (Cm) in the 
small parameter region between the synchronization and unsynchronization (the modular 
and global synchronization) grid points by varying Jinter or Moreover, to get more 

accurate transition curves, we divide a part of the parameter plane where the transition 
curves change rapidly into more minute grids and repeat the above computations. 

G. Sub- and Whole-population Statistical-Mechanical Spiking Measures 

We measure the degree of modular and global sparse synchronization in terms of realistic 
statistical-mechanical sub- and whole-population spiking measures, based on the ISPSR 
and the IWPSR kernel estimates and Rw{t) [7T] . Spike synchronization may be well 

visualized in the raster plots of spikes. For a synchronized case, the raster plot is composed of 
partially-occupied stripes (indicating sparse synchronization), and the corresponding ISPSR 
and IWPSR kernel estimates, Ri^\t) (/ = 1,2,3) and Rw{t), exhibit regular oscillations. 
Each Rh {i = 1, 2, 3,...) global cycles of Ri^\t) and Rw{t) begin from their left minimum, pass 
the central maximum, and end at the right minimum [also, corresponding to the beginning 
point of the next {i + l)th global cycles]; the 1st global cycles of Ri^\t) and Rw{t) appear 
after transient times of 10^ ms, respectively. Spikes which appear in the Rh global cycles 
of Ri^\t) and Rw{t) form the Rh stripes in the raster plots for the sub- and the whole- 
populations, respectively. To measure the degree of the sub- and the whole-population spike 
synchronization seen in the raster plots, statistical-mechanical sub- and whole-population 
measures and Ms^\ based on Ri^\t) and Rw{t), are introduced by considering the 
occupation pattern and the pacing pattern of spikes in the stripes for the sub- and the 
whole-populations, which corresponds to a simple extension of the case of non-modular 
networks PH. The sub- and the whole-population spiking measures (/ = 1, 2, 3) and 
of the Rh stripes [appearing in the Rh global cycles of Ri^\t) and Rw{t)] are dehned by 
the products of the sub- and the whole-population occupation degrees and of spikes 
(representing the density of the Rh stripes) and the sub- and the whole-population pacing 
degrees and P/’^^ of spikes (denoting the smearing of the Rh stripes), respectively: 

mP = OP ■ pP and Mp^ = Op^ ■ Pp\ (20) 
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The sub- and the whole-population occupation degrees and in the ith stripes are 
given by the fractions of spiking neurons in the ith stripes: 






and = 


N, 


is) 


M-L 


( 21 ) 


where and Npl are the numbers of spiking neurons in the zth stripes for the Jth 
sub-network and the whole network, respectively. For sparse synchronization with partially- 
occupied stripes, oP <C 1 and Op^ <C 1. The pacing degrees Pp'^ and Pp"' of sparse spikes 
in the ith stripes for the sub- and the whole-populations can be determined in a statistical- 
mechanical way by taking into account their contributions to the macroscopic ISPSR and 
IWPSR kernel estimates RP{t) and Rw{t), respectively. Instantaneous global phases ^P{t) 
of RP{t) and of Rw{t) are introduced via linear interpolation in the two successive 

subregions forming global cycles mi. The global phases ^P{t) and between the left 

minimum (corresponding to the beginning point of the Rh global cycle) and the central 
maximum are given by 
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and ^P{t) and $ 10 ( 2 ) between the central maximum and the right minimum [corresponding 
to the beginning point of the {i + l)th global cycle] are given by 


«>P(t) 


^w{t) 


27r{i — 1) -|- TT 


27r{i — 1) -|- vr 


t — t 


(max) 

lA 


Amin) Amax) 
^7,2+1 “ ^7,2 


t — t. 


(max) 


Amin) Amax) 

^W,i-\-l ^2R,2 


r_ Amax) Amin) 

lOI Ij i A ^ ^I,i+1 ’ 

r Amax) ^ . Amin) 

^w,i ^ ^w,i+l^ 


(24) 

(25) 


where and the beginning times of the Rh {i = 1, 2, 3, • • ■) global cycles of 

rPA) and Rw{t) [i.e., the times at which the left minima of RP{t) and Rw{t) appear in 
the Rh global cycles], respectively, and and tpA^ are the times at which the maxima 

of rPA) and Rw{t) appear in the Rh global cycles, respectively. Then, the contributions 
of the kth microscopic spikes in the Rh stripes occurring at the times 2j*], and tp\ to rP (t) 
and Rw{t) are given by cos<F^^^ and cos<F^“'\ where and are the global phases at 
the /cth spiking times [i.e., ^P = ^P{tpl) and = *h^(2^\.)]. Microscopic spikes make 
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the most constructive (in-phase) contributions to Ri^\t) and Rw{t) when the corresponding 
global phases and are 27rn (n = 0,1, 2,...), while they make the most destructive 
(anti-phase) contribution to Ri^\t) and Rw{t) when and are 27r(n — 1/2). By 
averaging the contributions of all microscopic spikes in the ith stripes to Ri^\t) and Rw{t), 
we obtain the pacing degrees and of spikes in the ith stripes: 
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(0 
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cos and ^ ’ 




(I) 
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S.. 


(26) 

k=l k=l 

where and are the total numbers of microscopic spikes in the zth stripes for the 


sub- and the whole-population, respectively. By averaging and of Eq. (20) over a 
sufficiently large number Ng of stripes, we obtain the sub- and whole-population statistical- 
mechanical spiking measures and Mg'^'^: 


Ns 


Ns 


= “““I = (27) 

^ i=l * i=l 

Here, we follow 3 x 10^ global cycles in each realization, and obtain average occupation 
degrees, average pacing degrees, and average statistical-mechanical spiking measures via 
average over 20 realizations. 


IV. EFFECT OF INTER-MODULAR CONNECTION ON FAST SPARSELY SYN¬ 
CHRONIZED RHYTHMS 


In this section, we investigate the effect of inter-modular connection on emergence of fast 
sparsely synchronized rhythms by varying both the inter-modular coupling strength Jinter 
and the average number of inter-modular links per interneuron in the clustered 

small-world network of inhibitory FS Izhikevich interneurons. In contrast to the case of 
non-modular networks, two kinds of sparsely synchronized states such as modular and global 
sparse synchronization are thus found. These sparsely synchronized states are characterized 


by employing diverse realistic measures, explained in Sec. Ill 


In the absence of inter-modular coupling, we consider three cases of the intra-modular 
dynamics in the Watts-Strogatz small-world sub-networks: (1) synchronized in all iden¬ 
tical sub-networks, (2) unsynchronized in all identical sub-networks, and (3) synchro¬ 
nized/unsynchronized in non-identical sub-networks. For each case, we study emergence 
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of sparsely synchronized popnlation states by changing both Jinter and for a hxed 

set of Inc, D, and (whose valnes are listed in Table I). In the snbsection 


WJi we start from the 1st case of intra-modular dynamics. To farther examine dependence 
on the type of intra-modular dynamics, we also study the 2nd and the 3rd cases in the 
subsection IIVBI 


A. 1st Case of Intra-Modular Dynamics: Synchronized in All Identical Sub¬ 
networks 


In this subsection, we consider the 1st case of intra-modular dynamics which are 
synchronized in 3 identical small-world sub-networks with the same rewiring probability 
Prewiring = 0.25. In the abseuce of inter-modular coupling (i.e., Jinter = 0), every sub¬ 
population in the small-world sub-networks exhibits identical sparse synchronization, as 
shown in Fig. Clear stripes are formed in the raster plot of Fig. [^a). The density of 

stripes is sparse because only a small fraction (about 0.22) of the total L (= 10^) neu¬ 
rons in the sub-population hre in each stripe. Due to presence of these sparse stripes, the 
ISPSR kernel estimate Ri^\t) shows fast regular oscillation with sub-population frequency 
~ 147 Hz, as shown in Figs. ^b)-^c). For the case of individual neurons, the inter¬ 
spike interval (ISI) histogram has multiple peaks appearing at multiples of the period Tj 
(~ 6.8 ms) of Ri^\t) (i.e., skipping of spikes occurs at random integer multiples of Tj) [see 
Fig. [^d)]. Because of this stochastic spike skipping (also called the stochastic phase locking) 
[I9H2I1 El Eol [ZH [83HH5], individual neurons exhibit stochastic and intermittent spike dis¬ 
charges, and hence partial occupation occurs in the stripes of the raster plot. In contrast to 
sub-population rhythms, the distribution of mean bring rates (MFRs) of individual neurons 
shows a peak near (~ 33 Hz) which is much less than the sub-population frequency 
[see Fig. |^e)]. In this way, bring activity of individual neurons dibers distinctly from the 
population oscillatory behavior for the case of sparse synchronization [alllMHlElllZQllII]. 
For more details on the sparse synchronization in the (non-modular) small-world network, 
refer to [64] . 

From now on, we employ the methods for characterizing population dynamics in Sec. in 
and investigate the ebect of inter-modular connection on sparse synchronization by chang¬ 
ing the inter-modular coupling strength Jinter for = 20. Figures [^a)-[^c) show the 
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raster plots of spikes in the three sub-populations for Jinter = 10, 500, and 2500, respectively. 
The corresponding ISPSR and IWPSR kernel estimates, and Rw{t) of Eqs. (10) and 

(12), for Jinter = 10, 500, and 2500 are also shown in Figs. |^dl)-|^f4), respectively. For 
small Jinter, the inter-modular coupling strength plays a constructive role to favor the pac¬ 
ing between spikes in each sub-network, as shown in the case of Jinter = 10. For each Jth 
sub-population, sparse stripes are formed in the raster plot and Ri^\t) shows a regular 
oscillation, as shown in Fig. [^a) and Figs. |^dl)-|^d4). Hence, each sub-population ex¬ 
hibits sparse synchronization. However, the intra-modular dynamics of sub-networks make 
some mismatching because both the stripes and the ISPSRs between the sub-networks are 
shifted. Vertical gray lines which pass minima of Rw{t) are drawn as reference lines for 
matching between Rlf\t) (/ = 1,2,3) [where the minima of Ri^\t) [Rw{t)] are denoted 
by solid (open) circles]. As a result of mismatching, the degree of sparse synchronization 
in the whole population becomes less than that in the sub-networks [i.e., the amplitude 
of Ru,{t) is less than that of Ri^\t)], and this kind of population behavior for Jinter = 10 
is referred to as the modular sparse synchronization because it reveals the modular struc¬ 
ture. With increasing Jinter, the mismatching degree between the intra-modular dynamics 
of sub-networks decreases, although the stripes in the raster plots become more sparse due 
to increased inhibition. Eventually when passing a threshold Jinter (— 268), intra-modular 
dynamics of sub-networks begin to make perfect matching. As a result, the population 
behavior becomes globally identical, independently of the cluster structure, as shown in 
Fig. Hb) and Figs. [^el)-[^e4) for Jinter = 500 [where all the minima of Ri^\t){I = 1,2,3) 
lie on the reference vertical line passing the minima of i?^„(t)], and it is referred to as the 
global sparse synchronization. However, for sufficiently large Jinter, due to strong inhibition 
the inter-modular coupling strength plays a destructive role to spoil the pacing between 
sparse spikes. Hence, as Jinter passes a higher critical value Jinter h (— 1657) the global 
sparse synchronization breaks into unsynchronization. As an example, refer to the case of 
Jinter = 2500. Sparse spikes in the raster plot of each sub-network are completely scattered 
without forming any stripes [see Fig. ^c)], and hence each ISPSR kernel estimate R^Pit) be¬ 
comes nearly stationary (i.e., every sub-network exhibits an unsynchronized state), as shown 
in Figs. ^fl)-^f4). We now vary not only Jinter but also and investigate emergence 

of modular and global sparse synchronization in the whole J inter -plane by using the 
method explained in the subsection HIF Thus, we obtain the state diagram in Fig. [^g). 
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Modular sparse synchronization emerges for small Jinter or in the “L”-shaped gray 

region, while in the dark gray region global sparse synchronization appears. For large Jinter 
(> 1572), unsynchronization occurs between the modular and the global synchronization. 
Changes in the population behaviors along the routes I, II, and III in Fig. |^g) are given in 
the following subsubsections. 


1. Effect of The Inter-Modular Coupling Strength on Population Synchronization along The 
Route I 


In order to study the effect of the inter-modular coupling Jinter on the population syn¬ 
chronization, we consider the case of the route I with = 20. Some results for this 

case are given for Jinter = 10, 500, and 2500 in Figs. |^a)-|^f4). As Jinter IS increased, a 
transition from modular sparse synchronization to global sparse synchronization when pass¬ 
ing a threshold Jinteri— 268), and eventually to unsynchronization when passing a higher 
threshold Jinterhi— 1657) occurs. The higher threshold Jfnterh f^r the transition to un¬ 
synchronization is determined through calculation of the sub- and the whole-population 
order parameters and {Pw)r of Eqs. (13) and (14), where {■ ■ ■)r denotes an average 


over realizations. Figures S a4) show plots of {O^P)r and {Oyf)r versus Jinter- For 
Jinter < Jinterh (— 1657), Synchronized states exist because the values of and {Oni)r 

become saturated to non-zero limit values for large L. When passing the higher threshold 
threshold Jfnterhi ^ transition to unsynchronization occurs because the order parameters 


and {Pw)r tend to zero as L —)■ cx). These unsynchronized states seem to appear due 
to a destructive effect of strong inhibition spoiling the pacing between sparse spikes. Here, 
we present two explicit examples for the synchronized and the unsynchronized states. First, 
we consider the synchronized case for Jinter = 1600. For L = 10^, sparse stripes are formed 
in the raster plot of spikes for each sub-network, and the ISPSR and the IWPSR kernel 
estimates R^P{f) and Rw{t) show regular oscillations, although there are some variations in 
the amplitudes [see Figs. |^b) and|^d)]. As L is increased to L = 10^, stripes in the raster 
plots become more clear, and RP\t) and Ruj{t) display more regular oscillations with nearly 
the same amplitudes, as shown in Figs. 0c) and0e ). Consequently, the population state for 
Jinter = 1600 seems to be synchronized because RP{t) and Rw{t) tend to show more regular 
oscillations as L goes to the inhnity. As a second example, we consider the unsynchronized 
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case of Jinter = 1700. As shown in Fig. [^f) for L = 10^, sparse spikes are scattered withont 
forming any stripes in the raster plot, and R^P{t) and Rw{t) in Fig. ^h) show little noisy 
flnctnations. In contrast to the synchronized case, as L is increased to L = 10^, sparse spikes 
become more scattered, and conseqnently RP{t) and Rw{t) become nearly stationary, as 
shown in Figs. |^g) and |^i). Hence, the popnlation state for Rnter = 1700 seems to be 
unsynchronized because RP{t) and Ru]{t) tend to be nearly stationary as L increases to the 
inhnity. 

In order to further understand the above synchronization-unsynchronization transition, 
we investigate the effect of inter-modular connection on the “microscopic” dynamical cross¬ 
correlations between neuronal pairs. As examples, we reconsider the same cases of Jinter = 
1600 and 1700 as in Fig. Figure [^al) shows the plots of the spatial cross-correlation 
functions of Eq. 17 versus I for L = 10^ in the case of Jinter = 1600. These spatial 
correlation functions are nearly non-zero constant (~ 0.04) in the whole range of /, 
and hence the correlation length r/j becomes L/2 (=500) covering the whole sub-networks 
(note that the maximal distance between neurons is L/2 because of the ring architecture 
on which neurons exist). Consequently, each sub-network is composed of just one single 
synchronized block. For L = 10^, the flatness of CP in Fig. ^a2) also extends to the whole 
range {I = L/2 = 5000) of the Jth sub-network, and hence the correlation length becomes 
rjj = 5000, which also covers the whole sub-network. In this way, for Jinter = 1600, due 
to a constructive role of Jinter favoring the pacing between sparse spikes, the correlation 
length rji seems to cover the whole sub-network, independently of L. For this case, the 
normalized correlation length fji (= ^), representing the ratio of the correlation length 
Tji to the sub-network size L (i.e., denoting the relative size of synchronized blocks when 
compared to the whole sub-network size), has a non-zero limit value, 1/2, and consequently 
synchronization emerges in each sub-network. In contrast, for Jinter = 1700 the spatial cross¬ 
correlation functions CP are nearly zero, independently of L, as shown in Figs. ibl)§:b2). 
For this case, due to a destructive role of Jinter spoiling the pacing between sparse spikes, 
the correlation length rji becomes nearly zero, and hence no synchronization occurs in each 
sub-network. 

We now investigate the type of synchronization through measurement of the matching 
degree between the intra-modular dynamics in sub-networks in the synchronized range of 0 < 
Jinter < ^interp— 1657) aloug the route I in Fig. |^g)- Figures |^al)-j^a5) show the raster 
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plots of spikes in the three sub-networks for Jinter =10, 30, 70, 400, and 1200, respectively. 
The ISPSR and the IWPSR kernel estimates, Ri^\t) and Rw{t), for Rnter =10, 30, 70, 400, 
and 1200 are also shown in Figs. |^bl)-|^b5), respectively. For each Jth sub-population, 
sparse stripes are formed in the raster plot of spikes and the ISPSR kernel estimate Ri^\t) 
shows a regular oscillation with global frequency fp^'^ ~ 147 Hz. Hence, each sub-population 
shows sparse synchronization. For the case of modular sparse synchronization for Jinter = 10, 
30, and 70, the intra-modular dynamics of sub-networks make some mismatching because 
both the stripes and the ISPSR kernel estimates between the sub-networks are shifted [see 
Figs. i:bl)-ib3) where the minima of R^f^ {t) (denoted by solid circles) lie off the reference 
vertical lines which pass the minima of Rw{t) (represented by open circles)]. Hence the 
amplitude of Rw{t) becomes less than that of Ri^\t). As Jinter is increased, the mismatching 
degree decreases, and hence the amplitude of Rwit) increases. Eventually, when passing a 
threshold J^er (— 268) global sparse synchronization occurs. Hence, for Jinter = 400 and 
1200, intra-modular dynamics of sub-networks (shown in their raster plots and ISPSR kernel 
estimates) make perfect matching [i.e., the minima of Ri^\t) he on the reference vertical 
lines, as shown in Figs. [^b4)-[^b5)], and hence the amplitude of Rw{t) becomes the same as 
that of Ri^\t). The matching degree between the intra-modular dynamics of subnetworks 
may be measured through calculation of the cross-correlation modularity measure Cm of 
Eq. 19 Figure [^d) shows the plot of {CM)r versus Jinter where {■ ■ ■)r denotes average 
over realizations. As Jinter is increased, {CM)r increases monotonically, and eventually when 
passing the threshold J*nter (— 268), its value becomes 1. Hence, for Jinter < Jinter modular 
sparse synchronization (with {CM)r < 1) emerges, while global sparse synchronization (with 
{CM)r = 1) appears for < J^er < J*nter,h- 

We also measure the degree of modular and global sparse synchronization in the synchro¬ 
nized range of 0 < Jinter < Jinter,h- showu in Figs.[^al)-|^a5), spike synchronization may 
be well visualized in the raster plots of spikes. For a synchronous case, “stripes” (composed 
of spikes and representing population synchronization) appear successively in the raster 
plot. For measurement of the degree the sub- and the whole-population spike synchro¬ 
nization seen in the raster plots, realistic statistical-mechanical sub- and whole-population 
measures and are introduced by considering the occupation pattern (represent¬ 

ing the density of the stripes) and the pacing pattern (denoting the smearing of the stripes) 
of spikes in the stripes for the sub- and the whole-populations, as explained in the sub- 
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section IIIG By varying Jinter) we follow 3 x 10^ stripes (i.e., 3 x 10^ global cycles) in 
each realization, and through an average over 20 realizations, we obtain the sub- and the 
whole-population occupation degrees {O^P)r and {Ow)r of Eq. (21), the sub- and the whole- 
population pacing degrees {Ps^'^)r and {Pw)r of Eq. (26), and the statistical-mechanical sub- 
and whole-population spiking measures and of Eq. (27), and the results are shown 
in Figs. [^al)-[^c4). For the case of modular synchronization [occurring on the left region 
of the vertical dotted threshold line for Jinter = Jinter (— 268)], both the occupation degree 
{Ow)r and the pacing degree {Pw)r for the whole-population are less than those for the sub¬ 
populations because of mismatching between the intra-modular dynamics of sub-networks. 
As Jinter E increased, their mismatching degrees become smaller, and eventually {Ou,)r and 
{Pw)r for the whole population become the same as those for the sub-populations for the 
case of global synchronization (occurring on the right region of the vertical dotted threshold 
line) due to perfect matching between the intra-modular dynamics of sub-networks. We 
hrst consider the occupation degree which characterizes the sparseness degree of population 
synchronization. For the sub-populations, the occupation degrees decrease mono- 

tonically because of increase in inhibition with increasing Jinter- In the case of modular 
synchronization, typical IWPSR kernel estimates Rw{t) show faster and smaller-amplitude 
oscillations with the whole-population frequency larger than the sub-population fre¬ 
quency fp\ and hence the occupation degree (0^)r for the whole-population becomes less 
than As Jinter E increased, {Ow)r increases and approaches {Oi^^)r due to decrease 

in the mismatching degree between the intra-modular dynamics of sub-networks, and even¬ 
tually when passing the threshold Jinter (i-®-; case of global synchronization) they 

become the same and then decrease with increasing Jinter- We note that modular and global 
synchronization is sparse one because is much less than unity [i.e., only a small 

fraction of the total L (= 10^) neurons in the sub-population Ere in each stripe]. Next, we 
consider the pacing degree between spikes in the stripes. For relatively small Jinter-, with 
increasing Jinter the sub-population pacing degree {Ps^'^)r increases due to a constructive 
role of Jinter favoriug the pacing between spikes, while for large Jinter {Ps^'^)r decreases as 
Jinter F increased because of a destructive role of Jinter spoiling the pacing between spikes. 
Through competition between these constructive and destructive roles of Jinter a “plateau” 
with high pacing degree is formed in a relatively wide region of intermediate Jinter for the 
case of global sparse synchronization. The whole-population pacing degree {Pw)r (which is 
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less than or equal to (Pi^^)^) also exhibits similar behavior. Consequently, both the sub- 
and the whole-population statistical-mechanical spiking measures and Ms^'^ (which are 
obtained by taking into consideration both the occupation and the pacing degrees of spikes 
in the stripes) show bell-shaped curves with their peaks at Jinter — 202 (corresponding to 
modular sparse synchronization) and 287 (corresponding to global sparse synchronization), 
respectively. For further understanding of the pacing degree between spikes, we also consider 
the spatial cross-correlations between neuronal pairs. Figures [^dl)-[^d5) show the spatial 
cross-correlation functions of Eq. (17) for Jinter =10, 30, 70, 400, and 1200, respectively. 
For the case of relatively small Jinter, with increasing Jinter the value of increases, but it 
decreases for large Jinter- For quantitative analysis, we introduce the sub-population spatial 
cross-correlation degree {{Cl^^)i)r given by double averaging of the spatial cross-correlation 
function over all lengths I and realizations. This sub-population spatial cross-correlation 
degree is a microscopic measure quantifying the cross-correlation degree between 

the microscopic IISR kernel estimates (t) without any explicit relation to the macroscopic 
occupation and pacing patterns of spikes. Figure [^e) shows plots of (obtained 

through average over 20 realizations) versus Jinter for / =1, 2, and 3. Similar to the case of 
the sub-population pacing degree {{Ci^^)i)r also display similar bell-shaped curves 

with peaks in the region of global synchronization. Hence, the statistical-mechanical pac¬ 
ing degree between spikes seems to be somewhat associated with the microscopic spatial 
cross-correlation degree between neuronal pairs. 


2. Effect of The Average Number of Inter-Modular Connections along the Routes II and III 

In addition to the above study on the effect of Jinter along the route I for = 20, we 

also investigate the effect of average number of inter-modular connections per interneuron 
on emergence of modular and global sparse synchronization along the routes II 
and III for Jinter = 500 and 2500, respectively [see Fig. [^g)]. For the case of the route II 
with Jinter = 500, the raster plots of spikes in the three sub-populations for = 2, 

5, 20, and 50 are shown in Figs. |^al)-[^a4), respectively. The corresponding ISPSR and 
IWPSR kernel estimates, and Ru,{t), for = 2, 5, 20, and 50 are also shown 

in Figs. |^bl)-|^b4), respectively. For each Jth sub-population, sparse stripes are formed 
in the raster plot and the ISPSR kernel estimate Ri^\f) shows a regular oscillation. As 
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is increased, more clear stripes appear in the raster plots of sub-networks, and hence 
the amplitudes of increase. Furthermore, with increasing the mismatching 

degree between the intra-modular dynamics of sub-networks decreases, and eventually when 
passing a threshold (~ 9) perfect matching occurs. Figure ^c) shows the plot 

of the cross-correlation modularity measure {CM)r of Eq. (19) versus Thus, for 


modular sparse synchronization with {CM)r < 1 emerges, while global 


sparse synchronization with {CM)r = 1 appears for . In this way, with 

increasing the pacing degree between spikes increases monotonically thanks to the 

increase in the degree of effectiveness of global communication between spikes. Hence, 
plays only a constructive role to favor the pacing between spikes in sub-networks as 
well as the matching between the intra-modular dynamics of the sub-networks, in contrast 
to dual roles of Jinter for the case of route I. Hence, unsynchronization does not appear. This 
constructive role of may be seen explicitly in Figs. ^dl)-^f2). We hrst consider the 

occupation degree which characterizes the sparseness degree of spike synchronization. For 
the case of modular synchronization, the sub-population occupation degree in the 

sub-networks decreases a little with increasing while remains nearly constant 

for the case of global synchronization. On the other hand, as is increased the whole- 

population occupation degree {Ow)r increases and approaches {Oi^'’)r because of decrease 
in the mismatching degree between the intra-modular dynamics of the sub-networks, and 
eventually when passing the threshold they become the same and then remain 

nearly constant with increasing Hence, the constant behavior of {Oi^^)r and {0^)^ 

for the case of global synchronization (which may occur because the average inhibition given 
to each neuron is the same for constant inter-modular coupling strength, independently of 
is in contrast to the monotonically-decreasing behavior of and {Ow)r for 

the case of route I [refer to Figs. 0 al)-S a4)]. Since {Oi^'^)r < 1 [i.e., only a small fraction 
of the total L (= 10^) neurons in the sub-population fire in each stripe], modular and 
global synchronization is sparse one. Next, we consider the pacing degree between spikes in 
the stripes. For both cases of modular and global sparse synchronization, with increasing 
both the sub- and the whole-population pacing degrees and {Pw)r increase 

monotonically due to a constructive role of favoring the pacing between the spikes, 

in contrast to the bell-shaped behavior for the case of the route I [refer to Figs, [^bl)- 
[^b4)]. Consequently, both the sub- and the whole-population statistical-mechanical spiking 
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measures Ms^'^ and (which are given by the products of the sub- and the whole- 

population occupation and pacing degrees) increase monotonically in both cases of modular 
and global sparse synchronization, which is also in contrast to the case of the route I [refer 
to Figs. 0:ci)^3 c4)]. To further understand the pacing degree between spikes in the stripes, 
we consider the sub-population spatial cross-correlation degree between neuronal 

pairs [given by double averaging of the spatial cross-correlation function of Eq. (17) over 
all lengths I and realizations]. Figure ^g) shows plots of (obtained via average 

over 20 realizations) versus for / =1, 2, and 3. Similar to the case of the sub¬ 
population pacing degree also displays monotonic increasing behavior. 

Hence, the statistical-mechanical pacing degree between spikes seems to be associated with 
the microscopic spatial cross-correlation degree between neuronal pairs, like the case of route 
I. 

We also investigate emergence of modular and global sparse synchronization by increas¬ 
ing along the route III for Jinter = 2500 (which is much larger than that for the 

case of route II). Unlike the case of the route II, for small a destructive effect to 

decrease the pacing degree between spikes occurs due to strong inhibition for Jinter = 2500, 
and hence when passing a lower threshold (~ 6) a transition from modular sparse 

synchronization to unsynchronization occurs. However, with further increase in a 

constructive effect of to favor the pacing between spikes and the matching between 

the intra-modular dynamics of sub-networks overcomes the destructive effect of strong in¬ 
hibition. Consequently, a transition to global sparse synchronization occurs when passing 
a higher threshold (~ 26). These results are well shown in Figs. |^al)-|^h). The 

raster plots of spikes in the three sub-populations for =1, 5, 20, 30, and 50 are shown 

in Figs. [^al)-|^a5), respectively. The corresponding ISPSR and IWPSR kernel estimates, 
and Rw{t), for =1, 5, 20, 30, and 50 are also shown in Figs. gbl)i^b5), 

respectively. For each Jth sub-population, sparse stripes are formed in the raster plot and 
R!k\t) shows a regular oscillation, except for the unsynchronized case of = 20 where 

spikes are scattered without forming stripes in the raster plot and R^P (f) is nearly stationary. 
Figures ici)§: c4) show the sub- and the whole-population order parameters oP and 
which determine a threshold for the synchronization-unsynchronization transition. In the 
region of < Mppp , both oP and tend to zero in the thermody¬ 

namic limit of L —)■ cxD, and hence unsynchronized states appear due to a destructive effect 








of strong inhibition. On the other hand, for or 


‘'syn,l ^ '‘'syn^h 

values of {OY^r and {OYjr become saturated to non-zero limit values for large L, and hence 
synchronized states exist. Particularly, for sparsely synchronized states 

appear due to a constructive effect of favoring the pacing between spikes. The type 

of sparse synchronization may be determined in terms of the cross-correlation modularity 
measure Cm of Eq. (19) which is shown in Fig. §:d) 


For modular 

sparse synchronization with (C'M)r < 1 (he., some mismatching between the intra-modular 
dynamics of sub-networks) emerges, while global sparse synchronization with {CM)r = 1 
(i.e., perfect matching between the intra-modular dynamics of sub-networks) appears for 
■ The degree of synchronization is also measured in terms of the occu¬ 
pation degrees, the pacing degrees, and the statistical-mechanical spiking measures in the 
sub- and the whole-populations, which are shown in Figs. Ii<ii)-§g2). We hrst consider 
the case of modular sparse synchronization. As Mgyn is increased, both the occupation 
degree {oY^r and the pacing degree {P^Yr in fhe sub-networks decrease due to a destruc¬ 
tive effect of strong inhibition for Jinter = 2500. In the whole-population, with increasing 
Mgyn the occupation degree {Ow)r increases and approaches {oY^r because of decrease 
in the mismatching degree between the intra-modular dynamics of sub-networks, and the 
pacing degree {Pw)r decreases like the case of {P^Yr- Thus, both the sub- and the whole- 
population statistical-mechanical spiking measures and (which are given by the 

products of the sub- and the whole-population occupation and pacing degrees) decrease as 
increases. On the other hand, for the case of global sparse synchronization which 
is similar to the case of route II, the constructive effect of favoring the pacing 

between spikes dominates. Consequently, with increasing both {Ps^Yr and {Pw)r 

increase monotonically, while both {oY^r and {Oyj)r remains nearly constant because the 
average inhibition given to each neuron is the same for constant inter-modular coupling 
strength Jinter, independently of ■ Consequently, as is increased, the sub- 

and the whole-population statistical-mechanical spiking measures and Ms^'^ show a 
monotonic increase. As in the case of route II, the modular and global synchronization is 
sparse one because {oY^r is much less than unity. Furthermore, the statistical-mechanical 
pacing degree between spikes in the sub-population seems to be associated with the spa¬ 
tial cross-correlation degree {{cYYi)r between neuronal pairs (obtained through average 20 
realizations), which is shown in Fig. |^h). For the case of modular sparse synchronization 
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(^i.e., )5 with increasing {{C\^’)i)r decreases monotonically due 

to a destructive role of strong inhibition, while for the case of global sparse synchronization 
(i.e., ), {{Ci^^)i)r exhibits a monotonic increase because of a constructive 

role of 


(inter)^ 


f (inter) 




B. 2nd and 3rd Cases of Intra-Modular Dynamics 


To further examine dependence of the inter-modular connection effect on the type of intra- 
modular dynamics, we consider the 2nd and the 3rd cases of the intra-modular dynamics: (2) 
unsynchronized in the absence of inter-modular coupling in all identical sub-networks with 
Prewiring = 0.05 and (3) uou-identical sub-networks where in the absence of inter-modular 
coupling, the 1st sub-network with Prewiring = 0.25 is synchronized, the 2nd sub-network 
with Prewiring = 0.15 is also Synchronized, but the 3rd sub-network with Prewiring = 0.05 
is unsynchronized. Figure lOJ^a) shows the state diagram in the Jinter-Ms^^^^^ plane for 


the 2nd case of intra-modular dynamics. This state diagram is similar to that for the 
1st case in Fig. [^g), except for the appearance of “L”-shaped region of unsynchronization 
for small Jinter or Beyond the unsynchronized region, modular and global sparse 

synchronization appears in the gray and the dark gray regions, respectively. For large 
Jinter (> 1402) unsyuchrouization also occurs for large Jinter (> 1402) between the modular 
and the global synchronization. We also make an intensive investigation of emergence of 
modular and global sparse synchronization by changing Jinter along the route of = 20 

(corresponding to the 1st route I for the 1st case). As in the 1st case, we obtain the raster 
plots of spikes in the three sub-populations and the corresponding ISPSR and IWPSR kernel 
estimates, and Rw{t) for representative values of Jinter =50, 200, 600, 1000, and 2000; 

for brevity associated figures are not presented. Unlike the 1st case 1 of the intra-modular 
dynamics, for small Jinter (=50) unsynchronization occurs because of suiaW prewiring (=0.05). 
For this case, the clustering coefficient is high, and hence partial stripes (indicating local 
clustering of spikes) seem to appear in the raster plots of spikes. Thus, the raster plots 
show zigzag patterns intermingled with partial stripes with diverse inclinations and widths, 
and hence spikes become difficult to keep pace with each other. Consequently, R^P{t) and 
Rwit) become nearly stationary. However, as Jinter is increased and passes a lower threshold 
Jinter,i (— ^^e inter-modular coupling strength Jinter plays a constructive role to favor 
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the pacing between spikes in each sub-network, and synchronized states with regularly- 
oscillating Ri^\t) and Rw{t) appear for Rnter = 200, 600, and 1000. On the other hand, 
for large Rnter, due to strong inhibition Rnter plays a destructive role to spoil the pacing 
between sparse spikes. Hence, when passing a higher threshold Jinterh (— 1402), a transition 
to unsynchronization occurs, as shown for Jinter = 2000. For this case, sparse spikes in the 
raster plots in each sub-network are scattered without forming stripes, and hence both (t) 
and Rw{t) become nearly stationary (i.e., every sub-network exhibits an unsynchronized 
state). Similar to the 1st case, one can consider additional routes with hxed values of Jinter 
(e.g., 750 and 2300). As shown in the state diagram of Fig. [Io|)a), unsynchronization occurs 
for small in contrast to the 1st case. However, when passing a threshold 

modular synchronization appears, and then the population behaviors are similar to those 
for the 1st case. 

Finally, we consider the 3rd case of non-identical sub-networks where in the absence 
of inter-modular coupling, the 1st sub-network with Prewiring = 0.25 is synchronized, the 
2nd sub-network with Prewiring = 0.15 is also synchronized, but the 3rd sub-network with 
Prewiring = 0.05 is unsyuchrouized. Thanks to a constructive role of Jinter favoring the 
pacing between spikes, a transition to synchronization occurs in the 3rd sub-network when 
passing a lower threshold Jinter i- With increasing the value of Jinter i decreases 

due to a constructive effect of to favor the pacing between spikes. For J < Jinterh 

the 3rd sub-network is still unsynchronized, while the 1st and the 2nd sub-systems are 
synchronized. As Jinter passes Jinterh Ord sub-network becomes synchronized, and then 
modular synchronization occurs due to mismatching between the intra-modular dynamics 
of sub-networks. Here, we consider the case of Jinter > J inter,i- Figure [l0|)b) shows the state 
diagram in the Jinter -plane for the case of Jinter > 1 and > 1 (i.e., the 

region where the 1st and the 2nd sub-networks are synchronized but the 3rd sub-network is 
unsynchronized is not shown). We note that this state diagram is nearly the same as that for 
the case 1 in Fig.j^g). Modular sparse synchronization occurs in the “L”-shaped gray region, 
while global synchronization appears in the dark gray region. Unsynchronization also occurs 
for large Jinter (> 1371) between modular and global synchronization. When compared 
with the 1st case of the intra-modular dynamics, the regions of modular synchronization 
and unsynchronization are a little enlarged, while the region of global synchronization is 
somewhat contracted. We make an intensive study on appearance of modular and global 
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synchronization by increasing Jinter from the threshold J*nterii— 0-2) for fhe 3rd snb-network 
along the ronte of = 20 (corresponding to the 1st ronte I for the 1st case). Similar 

to the 1st case, we obtain the raster plots of spikes in the three snb-popnlations and the 
corresponding ISPSR and IWPSR kernel estimates, Ri^\t) and Rw{t) for representative 
valnes of Rnter =10, 100, 500, 1000, and 2000; for brevity associated hgnres are not presented. 
For Jinter = 10, dne to a constrnctive role of Jinter favoring the pacing between spikes, 
sparse stripes are formed in each Jth snb-popnlation (J = 1,2,3). However, in contrast to 
the case 1 of identical intra-modnlar dynamics, the smearing degree of stripes is different, 
depending on the snb-popnlation. The stripes for the 1st snb-popnlation (with Prewiring = 
0.25) are relatively clear, while the stripes in the other 2nd and 3rd snb-popnlations (with 
Prewiring = 0.15 and 0.05, respectively) are more and more smeared. Hence, the amplitudes 
of the regularly-oscillating Ri^\t) decrease as I is increased. As Jinter is further increased 
(e.g., Ji„ter=100 and 500), the pacing degree of spikes in the stripes increases for each 
sub-population, although the stripes become more sparse. However, for large Jinter (e.g., 
Jinter = 1000), due to a destructive role of Jinter fhe pacing degree of spikes in the stripes 
begins to decrease. Eventually, when passing a higher threshold Jinter h (— 1372) a transition 
to unsynchronization occurs. Thus, for Jinter = 2000, spikes are scattered in the raster plot, 
and both R^P{t) and Rw{t) are nearly stationary. To study the effect of one may 

also consider another routes with hxed values of Jinter (e.g., 750 and 2300). With increasing 
along these routes, population behaviors are similar to those for the 1st case, as 
shown in the state diagram of Fig. [l0|(b). 

V. SUMMARY 

Sparsely inter-connected modular structures are found in both mammalian brain anatom¬ 
ical networks and human brain functional networks, as in other complex systems such as 
social, technological, and biological networks. Modular organization of the brain network 
shows the anatomical substrate for segregation of the brain into specialized sub-regions with 
particular functional tasks. These specialized informations of different features are also in¬ 
tegrated to produce new useful information as a whole. In this way, the brain network is 
organized via the interplay between segregation (specialization) and integration (binding). 
We note that these real brains, composed of sparsely linked clusters, are far more complex 
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than minimal non-clustered models snch as small-world and scale-free networks. To take 
into consideration the modnlar strnctnre of the real brain, we considered clnstered small- 
world networks of inhibitory FS internenrons, and investigated the effect of inter-modnlar 
connection on emergence of sparsely synchronized rhythms at the snb- and whole-population 
levels by employing diverse realistic measures. By changing both the inter-modular coupling 
strength Jinter and the average number of inter-modular links per interneuron we 

made intensive study on emergence of sparsely synchronized population states along the 
three routes in the Jinter ~plane for the 1st case of the intra-modular dynamics: (1) 
synchronized in all identical sub-networks. Consequently, two kinds of sparse synchroniza¬ 
tion such as modular and global synchronization have been found to appear, in contrast 
to the case of non-clustered networks. Our main hndings are that the type and degree of 
sparse synchronization depend on the inter-modular parameters, Jinter and For 

lower values Jinter acts to favor the pacing between spikes, while for higher values it tends to 
spoil the pacing between spikes due to strong inhibition. On the other hand, with increas¬ 
ing it acts to monotonically increase the pacing between spikes, which results from 

increase in the degree of effectiveness of global communication between spikes. To examine 
dependence on the intra-dynamics in the sub-networks, we have also considered two other 
cases for the intra-dynamics: (2) unsynchronized in all identical sub-networks and (3) syn¬ 
chronized/unsynchronized in non-identical sub-networks. Figures [^g),[T^a) and[IO|(b) show 
the state diagrams, representing main features on the population states, for the 1st, 2nd, 
and 3rd cases, respectively. For the 2nd case, an “L”-shaped region of unsynchronization 
appears for small Jinter or because the intra-modular dynamics in all sub-networks 

are unsynchronized in the absence of inter-modular coupling. Beyond this “L”-shaped re¬ 
gion, the structure of the state diagram is similar to that for the 1st case. In the 3rd case, 
its state diagram is nearly the same as that for the 1st case. Due to non-identicalness of 
sub-networks, modular synchronization persists a little more, and hence the region of mod¬ 
ular synchronization is a little enlarged when compared to the 1st case. Moreover, in the 
case of modular synchronization the pacing degree of spikes varies depending on the sub¬ 
networks, in contrast to the 1st case of identical sub-networks. From the results in these 
three cases, it follows that the effect of inter-modular connections seems to be essentially 
the same, independently of the type of intra-modular dynamics in sub-networks. Finally, 
since changes in the coupling strengths and the synaptic connections are closely interwoven 
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with the brain plasticity [86], we expect that our results on the inter-modular connection 
effect in modular networks might have implications for the role of the brain plasticity in 
some functional behaviors related to population synchronization. However, explicit study 
on the inter-relation between inter-modular connection, population synchronization, and 
brain plasticity is beyond our present subject, and it is left as a future work. 
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TABLE I: Parameter values used in our computations for Figs. igg units of the capacitance, the 
potential, the current, and the time are pF, mV, pA, and ms, respectively. 

(1) Single Izhikevich FS Interneurons |74) 


C = 20 Vr = -55 vt = -40 Vp = 25 Vb = -55 

k = 1 a = 0.2 b = 0.025 c = -45 d = 0 


( 2 ) 


(3) 


(4) 


(5) 


External Stimulus to Izhikevich FS Interneurons 
Idc = 1500 D = 500 

Inhibitory GABAergic Synapse m 

Tl = I Tr = 0.5 Td = 5 Vsyn = -80 

Intra-modular Coupling in Small-world Sub-networks 

= 50 Jintra = 1400 Prewiring ■ Varying 

Inter-modular Connection 

: Varying Jinter ■ Varying 
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FIG. 1; Schematic representation of the clustered small-world network topology. The whole 
network consists of M(= 3) Watts-Strogatz small-world sub-networks, each of them containing 
L(= 20) interneurons. Within each sub-network, the average number of intra-modular synaptic 
inputs per interneuron is 4), while there are 8 sparse random inter-modular connections 

(denoted by heavy lines) between the sub-networks. 
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FIG. 2: Fast sparse synchronization in identical small-world sub-networks with Prewiring = 0.25 
(case 1) in the absence of inter-modular coupling, (a) Raster plot of neural spikes and (b) in¬ 
stantaneous sub-population spike rate kernel estimate Ri^\t). (c) One-sided power spectrum of 
ARi^\t) [= Ri^\t) — with the mean-squared amplitude normalization, (d) Inter-spike 

interval (ISI) histogram for the individual interneurons; the vertical dotted lines denote integer 
multiples of the global period Tj (~ 6.8 ms) of Ri^\t). (e) MFR (mean firing rate) distribution of 
individual interneurons. 
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FIG. 3; Modular and global sparse synchronization for the 1st case of intra-modular dynamics 
with Prewiring = 0.25. = 20 in (a)-(f4). Raster plots of neural spikes in the sub-networks 

(/ = 1,2,3) for Jinter = (a)10, (b) 500, and (c) 2500. Instantaneous sub- and whole-population 
spike rate kernel estimates, (I = 1, 2,3) and Rw{t)- Jinter = (dl)-(d4) 10, (el)-(e4) 500, and 

(fl)-(f4) 2500. Vertical gray lines for Jinter = 10 and 500 pass minima of Rw{t), and the minima of 
Ri^\t) [R«,(t)] are denoted by solid (open) circles, (g) State diagram in the Jinter- plane. 
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FIG. 4: Realistic thermodynamic order parameters for measurement of the threshold for the 
synchronization-unsynchronization transition along the route I with = 20 for the 1st 

case of intra-modular dynamics with Prewiring = 0.25. (al)-(a4) Plots of the sub- and the whole- 
population order parameters logj^Q (I = 1,2,3) and logiQOw versus Jinter- Sparse synchro¬ 
nization for Jinter = 1600: raster plots of neural spikes in the three sub-networks for (b) L = 10^ 
and (c) L = 10^ and instantaneous sub- and whole-population spike rate kernel estimates 
(/ = 1,2,3) and Rw{t) for (d) L = 10^ and (e) L = 10^. Unsynchronization for Jinter = 1700: 
raster plots of neural spikes in the three sub-networks for (f) L = 10^ and (g) L = 10^ and instan¬ 
taneous sub- and whole-population spike rate kernel estimates (/ = 1,2,3) and Rw{t) for 

(h) L = 10^ and (i) L = 10^. 
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FIG. 5; Spatial cross-correlation functions for the synchronized and the unsynchronized states along 
the route I with = 20 for the 1st case of intra-modular dynamics with Prewiring = 0.25. 

Sparse synchronization for Jinter = 1600 with the spatial correlation length rjj covering the whole 
sub-network; spatial cross-correlation functions in the three sub-networks (/ = 1,2,3) for 
(al) L = 10^ and (a2) L = 10^. Unsynchronization for Jinter = 1700 with r// ~ 0; spatial 
cross-correlation functions in the three sub-networks for (bl) L = 10^ and (b2) L = 10^. 
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FIG. 6: Realistic cross-correlation modularity measure for determining the threshold for modular- 
global sparse synchronization transition along the route I with = 20 for the 1st case of 

intra-modular dynamics with Prewiring = 0.25. Raster plots of neural spikes in the sub-networks 
(/ = 1,2,3) for Jinter = (al) 10, (a2) 30, (a3) 70, (a4) 400 and (a5) 1200. Instantaneous sub- and 
whole-population spike rate kernel estimates (/ = 1,2,3) and Rw{t) for Jinter = (bl) 10, (b2) 
30, (b3) 70, (b4) 400 and (b5) 1200. Vertical gray lines pass minima of Rw{t), and the minima of 
Ri^\t) [Rw{t)] are represented by solid (open) circles. Temporal cross-correlation functions 
between the instantaneous sub-population spike rate kernel estimates Ri^\t) and Ri'^\t) of the 
sub-networks I and J for Jinter = (cl) 10, (c2) 30, (c3) 70, (c4) 400 and (c5) 1200. (d) Plot of the 
cross-correlation modularity measure {CM)r versus Jinter- 
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FIG. 7: Realistic statistical-mechanical spiking measure for measurement of the degree of modular 
and global sparse synchronization along the route I with = 20 for the 1st case of intra- 

modular dynamics with Prewiring = 0.25. Vertical dotted lines in (al)-(c4) and (e) represent the 
threshold of J = Jitter (— 268). (al)-(a4) Plots of the sub- and the whole-population occupation 
degrees and {0^)^ versus Jinter- (bl)-(b4) Plots of the sub- and the whole-population 

pacing degrees and {Pw)r versus Jinter- (cl)-(c4) Plots of the sub- and the whole-population 

statistical-mechanical spiking measures and versus Jinter- Spatial cross-correlation 

functions in the three sub-networks for Jinter = (dl) 10, (d2) 30, (d3) 70, (d4) 400 and (d5) 
1200. (e) Plot of the spatial cross-correlation degree versus Jinter- 
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FIG. 8: Modular and global sparse synchronization along the route II with Jinter = 500 for the 
1st case of intra-modular dynamics with Prewiring = 0.25. Raster plots of neural spikes in the 
sub-networks (/ = 1,2,3) for = (al), (a2) 5, (a3) 20 and (a4) 50. Instantaneous sub- and 

whole-population spike rate kernel estimates Ri^\t) (/ = 1,2,3) and Rw{t) for = (al), 

(a2) 5, (a3) 20 and (a4) 50. Vertical gray lines pass minima of Rw{t), and the minima of Ri^\t) 
[R^(t)] are represented by solid (open) circles, (c) Plot of the cross-correlation modularity measure 
{CM)r versus Vertical dotted lines in (dl)-(g) denote the threshold (~ 9). 

(dl)-(d2) Plots of the sub- and whole-population occupation degrees and {0^)^ versus 

(el)-(e2) Plots of the sub- and the whole-population pacing degrees and {Pw)r 

versus (fl)-(f2) Plots of the sub- and the whole-population statistical-mechanical spiking 

measures and {Mg^'^)^ versus (g) Plot of the spatial cross-correlation degree 

{{CP)i)r versus 
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FIG. 9: Modular and global sparse synchronization along the route III with Jinter = 2500 for the 
1st case of intra-modular dynamics with prewiring = 0.25. Raster plots of neural spikes in the 
sub-networks (I = 1,2,3) for = (al) 1, (a2) 5, (a3) 20, (a4) 30 and (a5) 50. Instantaneous 

sub- and whole-population spike rate kernel estimates (I = 1, 2, 3) and Rw{t) for = 

(al) 1, (a2) 5, (a3) 20, (a4) 30 and (a5) 50. Vertical gray lines pass minima of Rwit), and the 
minima of Ri^\t) [Rw{t)] are represented by solid (open) circles. (cl)-(c4) Plots of the sub- and 
the whole-population order parameters log;^g {O^P)r {I = 1,2,3) and log;^o {^w)r versus 
In (d)-(h), break symbols are given in the unsynchronization regions; the left (right) parts of the 
break symbols correspond to the regions of modular (global) synchronization, (d) Plot of the 
cross-correlation modularity measure (C'M)r versus (el)-(e2) Plots of the sub- and the 

whole-population occupation degrees and (0«,)^ versus (fl)-(f2) Plots of the sub- 

and the whole-population pacing degrees {Ps^'^)^ and {Pw)r versus (gl)-(g2) Plots of the 

sub- and the whole-population statistical-mechanical spiking measures (Afi^^)^ and versus 

(h) Plot of the spatial cross-correlation degree versus Mg^^^K 
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FIG. 10: (a) State diagram in the Jinter-Mg^'^^'^ plane for the 2nd case of the intra-modular 
dynamics (i.e., in the absence of inter-modular coupling, identically unsynchronized sub-networks 
with Prewiring = 0.05). (b) State diagram in the plane for the 3rd case of the 

intra-modular dynamics (i.e., in the absence of inter-modular coupling, synchronized in the 1st 
sub-network with Prewiring = 0.25, less synchronized in the 2nd sub-network with Prewiring = 0.15, 
and unsynchronized in the 3rd sub-network Prewiring = 0.05). 
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